Exchange— correlation effects in the monocUnic to tetragonal phase 
stabilization of Yttrium— doped Zr02- a first— principles approach 

Davide Sangalli^'^ and Alberto Debernardi^ 

^MDM Lab - IMM - CNR via C. Olivetti, 
2 1-20864 Agrate Brianza (MB) Italy, European Union 
^ Consorzio Nazionale Interuniversitario per le Scienze dei Materiali ( CNISM) 

(Dated: December 16, 2011) 

Abstract 

We describe, within an ab-initio approach, the stabihzation of the tetragonal phase vs. the 
monoclinic one in Yttrium-doped Zirconia. The process is believed to be influenced from different 
mechanisms. Indeed we show that there is a delicate balance between the change in electrostatic 
and kinetic energy and exchange-correlation effects. In the tetragonal phase the perturbation 
induced by doping is better screened at the price of sacrificing correlation energy. Our work opens 
the opportunity to use the same approach to predict the tetragonal phase stabilization of materials 
like Zirconia or Hafnia, with different and less characterized dopants. 

PACS numbers: 64.70.K-, 71.15.Mb, 81.05.Je, 81.30.Bx 
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I. INTRODUCTION 



Zirconia {Zr02) is a hard and usually colorless material with a wide range of technological 
applicationsi. Being corrosion resistant it is used as a dental material and because of its 
low cost, durability, and close visual likeness to diamond, it is widely used to synthesize 
artificial gems. It is used as a thermal barrier in coating engines due to its high resistance. 
The addition of cations (as for example Y^~^) induces the generation of oxygen vacancies for 
charge compensation which makes it useful as oxygen sensor. Moreover Zirconia is an high 
dielectric constant (high-/?) material with potential applications in the micro-electronics. 
Finally, very recently, transition-metals doped Zirconia has been predicted to be a dilute 
magnetic semiconductor (DMS) with high Curie Temperature Tc with potential applications 
in the field of spintronics^; while pure and doped Zr02 has been proposed as a candidate 
material for resistive switching memories devices^"^ (ReRAM) exploiting the high vacancies 
mobility of the system. 

Pure Zr02 exhibits at ambient pressure three polymorphisms. The monoclinic (M) phase 
is stable at low temperature and is the less symmetric structure with the Zr^~^ ions exhibiting 
sevenfold coordination. Between 1440 and 2640 K the tetragonal (T) phase, with eightfold 
coordinated Zr^~^ ions is stable. Finally above 2640 K, till melting temperature (~ 3000K), 
the most symmetric cubic (C) phase is stabilized^*^. The only difference between the (C) 
and the (T) phase is a distortion of the oxygen sub-lattice with a spontaneous symmetry 
breaking. The (T) and (C) phases are more used in technological applications and both the 
(T)/(C) structures and the (T) — )• (C) doping induced phase transition (DIPT) have been 
well characterized in the literature both form the theoretical and the experimental point of 
view^i^"— , with Y^~^ ions the most used and studied dopants. 

On the other hand the (M) — )■ (T) DIPT is much less characterized, especially from the 
theoretical point of view. It presents a volume change of about 3 — 4% that causes exten- 
sive cracking in the material. This behavior destroys the mechanical properties of fabricated 
components and makes pure Zirconia useless for structural or mechanical application. More- 
over the (T) phase, is metastable in pure and lightly doped Zr02 over a very long time; a 
growth sample of doped Zirconia must be annealed in order to check if the reached tetrag- 
onal phase is stable or meta-stable. Hence a better understanding of the (M) — )■ (T) DIPT 
is desirable. In this paper we address the problem from a first principle perspective and, 
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Table I. Theoretical and experimental cell parameters of 
Zr02 in the monoclinic and in the tetragonal phase. 



Present work Theory Experiment 
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1.012 


1.012 


c/a 
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1.023 


1.032 
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99° 10' 


99°5' 


99° 14' 


99° 14' 






(T) Zr02 






a (A) 


5.11 


5.02 5.03 




5.07 


c/a 


1.030 


1.014 1.017 


1.026 


1.018 



in particular, we show that the DIPT is a balance of the mean field (MF) description (i.e. 
kinetic plus electrostatic: H = T + y^^^^ _)_ Vh'[p]) with exchange-correlation effects which 
cannot be captured by simplified models. 

II. FIRST PRINCIPLES DESCRIPTION. 
A. Computational details. 

We work within density functional theory^i^S (DFT) in the generalized gradient approx- 
imation^ (GGA) with ultrasoft pseudopotentials^"^. Both for Yttrium an Zirconium the 
pseudopotentials includes semicore electrons and non-local core correction. The Quantum 
Espresso^^ package is used to solve the Kohn-Sharr>22 (KS) equations for a super-cell with 
96 atoms (down to 92 when oxygen vacancies are considered). We checked that a cut-off 
of 35 Ry for the wave-functions and 400 Ry for the augmentation charge, a k-point grid 
2x2x2 are needed in order to have the energy difference between the tetragonal and the 
monoclinic phase converged to up to 5.0 * 10~^ eV per ZrO^ molecular unit (m.u.), where 
the energy difference is of the order of 0.1 - 0.01 eV/m.u. . 

As starting point we relaxed atomic positions and structure of pure Zr02 with 12 atoms 
in the unit cell of the (T) and (M) phases. The (C) phase is unstable and is not considered 
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Table II. Theoretical and experimental internal structural parame- 
ters of Zr02 in the WyckofS^ notation. For the tetragonal struc- 
ture the only free parameter is the z coordinate of the oxygen atoms 



(0.0, 


0.5, 0.25-4) 






Present work 


Theory 


Experiment 






(M) Zr02 








Ref.fS] 


Ref.[14] 


Zr 


(0.276,0.044,0.210) 


(0.277,0.043,0.210) 


(0.275,0.040,0.208) 


Oi 


(0.065,0.327,0.350) 


(0.064,0.324,0.352) 


(0.070,0.332,0.345) 


On 


(0.451,0.757,0.475) 


(0.450,0.756,0.479) 


(0.450,0.757,0.479) 






(T) Zr02 


Ref.[13] Ref.[26j 


dz 


0.057 


0.044 


0.057 0.042 



at zero doping^. The results, reported in Tables III Atlll Al are in good agreement with 
experimental data and previous theoretical works within the DFT error (?5i 1 — 2%). 

We correctly find that the (M) phase is favored with an energy difference AEm~t = 
0.109 [eV/mol] in agreement with previous works (0.063 [eV/mol]^, 0.144 [eV/mol]^); the 
experimental estimation is 0.063 [eV/mol}^. 



B. Yttria stabilized Zirconia 

The (T) and the (C) phase of Zirconia can be stabilized with Yttria (Y2O3) doping, with a 
phase diagram for Yttria-stabilized Zirconia (YSZ) which has been extensively characterized 
experimentally. Also the (T) and (C) phases of Y-doped Zirconia have been studied within 
a theoretical approach in a number of works^"— , while no systematic study of the Y-doped 
(M) phase or of the (M) — > (T) DIPT has been reported in the literature. 

F-doping is known to induce oxygen vacancies. We show in Fig. [T] how the Kohn-Sham 
(KS) density of states changes if vacancies are present both in pure and doped Zr02 crystals 
at low doping^^. 

Moreover we compute the energy the system gains producing vacancies (see the caption 
of Fig. [T]for more details) for different doping concentrations. A correct modelling of the ma- 

4 



6 7 8 9 10 11 12 13 14 15 16 17 18 
energy [eV] 



6 7 8 9 10 11 12 13 14 15 16 17 



energy [eV] 



(a)DOS Zr02 with Y doping 



(b)DOS Zr02 with vacancy 





























•1 ^ 'V \ 
'/ ' 1 1 






A 




i 








fit ' 
/'■ 
/' ■' 


^' iVi 


'i 

'// 

V 


],],],], 






/' ' 
// / 


1 ,v 



Y doping 



9 10 11 12 13 14 15 16 17 18 
energy [eV] 




(c)DOS Zr02 with Y2O3 doping 



(d) Energy gain per vacancy 



Figure 1. (colors online) Kohn-Sham density of states for (a) Yttrium doped zirconia (Zri_a;02l^), 
(6) zirconia with oxygen vacancies {Zr02-yVy) and (c) zirconia with one vacancy each two Yttrium 
atoms {Zri-x02-yYxVy) in the tetragonal phase. The DOS are at x = 0.25 and y = 0.125 (red 
continuous line) and x = 0.125 and y = 0.0625 (orange dashed line); the DOS of pure Zr02 
(thinner black dot-dashed line) is also shown. The system doped with Yttrium gains energy when 
vacancies are created, releasing oxygen molecules, (panel (d)). Here Egain = ^ / y{E{Zri-x02Yx) — 
[E{Zri_x02-yYxVy) + y/2/i(02)]}, with x = 2y for the "Y-doping" case (on the left) and x = 
for pure Zr02 (on the right). The creation of vacancies is instead highly unfavored in the undoped 
system. 



terial must consider the possible relative positions of the dopants (as substitutional defects) 
and the vacancies. We found out that the relative position of Y atoms among themselves 
plays a minor role, in agreement with the other works^, and in our system we kept these as 
far as possible to mimic uniform doping. Instead the position of vacancies respect to Y atoms 
or to other vacancies influences significantly the total energy of the system. Stapper et al.- 
reported that in the cubic system isolated vacancies tend to be next nearest neighbor (NNN) 
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to Y atoms with an energy gain of Ai^c ~ 0.34 eV per vacancy against the nearest neighbor 
(NN) configuration in a supercell with 1 Y atom and 1 vacancy. This resuh is confirmed also 
by Ostanin et al.^ in a supercell with 2 Y atoms and 1 vacancy with AE^. ^ 0.3 eV between 
the two configurations (O-Yl, O-Y2) = {NN, NNNN) and (O-Fi, O-Y2) = {NNN, NNN) 
in favor of the second in the (C) phase and AEf ^ 0.2 eV in the (T) phase. Also in the 
present work we found AEt ~ 0.12 eV considering 1 vacancy and 2 Y atoms. 

In the (M) phase two nonequivalent oxygen atoms exist, one with coordination 3 and an- 
other with coordination 4, and, accordingly there are four possible configurations, we found 
that the case (A^A^, A^A^A^A^)x3 favored of 0.21, 0.36 and 0.80 eV against (A^A^A^, A^A^A^)x3, 
{NNN, NNN)x4 and (A^A^, A^A^A^A^)a;4 respectively; here the labels x3 and x4 refer to the 
coordination that would have the vacant oxygen. At the best of our knowledge, no results 
have been reported in the literature about the position of vacancies in the (M) phase. 

The most favored configuration, for both the (T) and the (M) phase, is represented in 
Fig. El together with the forces acting on each atom in the initial geometry of pure Zr02- 
The Y ions are not attracted by the oxygen vacancies and indeed in the (M) phase they 
are pushed away from the vacancy. Hence both the analysis of the energetically favored 
configuration and the internal relaxation of the ions does not follow the intuitive scheme of 
charged defects with the vacancies a +2e and the Y ion a —1 charged sites, suggesting that 
exchange-correlation effects play a role in the stabilization mechanism. 

Also the relative vacancy-vacancy orientation and position is in contrast with an elec- 
trostatic based model. In the (C) / (T) phase, while Stapper et al.- based their calculation 
on the assumption that oxygen vacancies should remain as far as possible, Ostanin et al.- 
results suggest that oxygen vacancies tend to couple along the [111] direction with a single 
Zr ion in between, thus at a relative distance of about ~ 4.47 A. Also experimental X-ray 
datar^ seems to support this idea. No data about the orientation of oxygen vacancies in the 
(M) phase is present in the literature. 

Changing the relative position of di-vacancy complex in our model (i.e. two vacancy at 
X = 12.50%) we found differences in total energy of AEtot ~ 1.1 eV (per di-vacancy). The 
oscillations in the mean field (MF), i.e. the sum of the kinetio^^ and electrostatic energy, 
and the exchange-correlation (xc) component of the energy are much more pronounced, 
with AExc ~ 8.2 eV and AEmf ~ 7.2 eV. We will discuss this more in detail in the next 
section. Our results are, in some aspects, different from the ones suggested by the two cited 
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(a) Tetragonal (b) Monoclinic 



Figure 2. (colors online) Tetragonal (a) and monoclinic (b) structure of Zr(^i_2x)Y2xO(^2-x)Vx, most 
stable configuration at x = 0.03125. The forces acting on the ions in the "starting geometry", i.e. 
with the atoms placed at the coordinates of the relaxed Zr02 system, are represented with blue 
arrows. In the tetragonal structure the oxygen atoms (O, small red spheres) closer to the oxygen 
vacancy (V represented as a small black sphere in the figure) move toward it of ~ 0.3 — 0.5 A, 
while the Zirconium atoms {Zr, big dark/magenta sphere) next nearest neighbor move away from 
it of ~ 0.1 — 0.2 A from the "starting geometry". The Yttrium atom {Y, biggest light /yellow 
sphere) is almost fix. The behavior is similar in the monoclinic structure, though the vacancy has 
just 3 nearest neighbor here, one of which is an Yttrium atom. The Yttrium atom moves away 
from the vacancy suggesting that a charged-defects interaction model based on electrostatic is not 
correct and that exchange-correlation effects play a major role. 

works. Indeed we found that, while vacancies does not repel, the configuration suggested as 
most favorable for the (T) structure by Ostanin et al.^, is not the lowest in energy. Indeed, 
surprisingly, the latter comes out to be unstable, with one of the two vacancies which changes 
position if we let the system relax (Fig. [3]). This last result suggests that, in presence of 
vacancies, oxygen atoms can move with (nearly) no potential barrier to overcome. This is a 
theoretical evidence that YSZ is a good ionic conductor. Indeed this property is important 
for many applications and has been investigated in other works^-. 

We also considered some cases at atomic doping concentrations of x = 18.75% and one 
test case at x = 25.00%. Both for x = 12.50% and x = 18.75% there is a huge number of 
possible relative position of the oxygens and we did not try to systematically explore which 
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(a)Starting geometry (b)Relaxed geometry 



Figure 3. (colors online) Starting (a) and final relaxed geometry (6) of the Zri^x02-x/2^xyx/2 
system with x = 12.5%. O atoms are dark (red) small spheres, Zr atoms are big dark spheres 
(magenta), and Y atoms are biggest light spheres (Yellow). The oxygen which changes position 
is a small light (light red) sphere. Before relaxation the two vacancies are aligned along the [111] 
direction at a distance of 8.45 Bohr; they can be detected as in the figure the oxygens are aligned in 
groups of four atoms everywhere except where a vacancy exist. After the relaxation two vacancies 
are aligned along the [101] direction at a distance of ~ 6.2 Bohr (this is the distance that the 
missing oxygens would have at the geometry of pure Zr02)- 



is the best configuration. At a; = 12.25% the best found configuration has two vacancies 
aligned in the xy plane. A simple explanation of this fact could be provided considering the 
axial anisotropy of the dielectric constant (k) of the (T) phase. The component hi^x = i^yy 
are larger than the Kzz, thus providing a more effective screening of charged oxygen vacancies 
placed along the xy plane. However we did not find this to be the sole mechanism, as some 
of the other configurations with the vacancies aligned in the xy plane have higher energy 
than some configurations with the vacancies aligned, for example, in the xz plane. 

To conclude this section, the behavior of the oxygen vacancies present some differences, 
and some similarities, between the two phases of Zirconia. A key important difference 
however is that the formation of vacancies is more favored in the (T) phase than in the (M) 



phase as one can see from Fig. 1(d) This fact indeed implies that 1^203 doping tends to 
stabilize the (T) phase. 
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Figure 4. (colors online) The difference in total energy per molecular unit (up panel) between 
the monoclinic and the tetragonal phase decrease increasing the Y2O3 doping (black circles). The 
result is not simply the sum of the effect induced by Y doping without vacancies (violet triangles, 
dot dashed line for the linear fit) and oxygen vacancies without Yttrium (green squares, dashed 
line for the linear fit). When no Yttrium is present the value of the x coordinate is chosen so 
that the number of vacancies for the Y2O3 doping case is equal to the number of vacancies for the 
oxygen vacancies only case. In the bottom panels the variations in the total energy as a function 
of the configurations is compared with the variations in the mean field and exchange-correlation 
energy. The pink shadowed area are a guide for the eyes, to highlight the trend of the two energy 
components. 



C. Theoretical stabilization: xc— energy and vacancies. 



The configuration explored allow to model the (M)— )'(T) DIPT as a function of the Yttria 
doping. In Fig. H] we see that the energy difference between the (T) and the (M) phase 
decrease increasing the doping concentration. This is a clear signature that we are correctly 
describing the DIPT. However the stabilization of the (T) phase happens, after a linear fit 
of our data, at higher doping concentration, ^ 14%, than the experimental value ~ 7%r^. 
The main error is due to an overestimation of the zero-doping energy difference of 0.056 
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eV/m.u. {Eexp ~ 0.063 eV/m.u.— while Edft ~ 0.109 eV/mu.u.). Indeed it is reasonable 
to assume that the trend of the energy difference is better computed than its absolute value 
and accordingly assuming a constant "zero-doping error" for every Y concentration we can 
subtract it (the zero energy in Fig. H] shifts from the black line to the blue line) obtaining 
xdipt ~ 7.5%, which is in excellent agreement with the experimental value. Similar results 
are obtained if one consider, instead of the linear fit, the energy difference between the two 
best configurations (orange dashed line in Fig. H]). 

For the pure Zr02 we have also computed the contribution to the Helmholtz free energy 
due to thermal excitations of phonons. We found that the difference in the thermal con- 
tribution is negligible at room temperature, while the difference of the zero point phonon 
energy between the two phases is ~ 0.005 eV/m.u. , with a small effect on the value of xbipt 
(In Fig. mthis amounts to the shift of the zero from the dashed to the continuous line). 

It is interesting to compare how the mean field energy, Emf[p\ = T[p] + V^^^ + Vh[p], 
and the exchange-correlation energy, E^cip]-, change considering different configurations. In 
Fig. m bottom panels, we see that the two component oscillate about 10 times more than 
the total energy as a function of the chosen atomic configuration. The mean field energy, 
in pure Zr02, is lower in the (T) phase (Fig. HJ central bottom panel); that is, if the xc- 
energy is neglected, the (T) phase is thermodynamically the most stable. The inclusion of 
Y doping and oxygen vacancies enhance this aspect as, on average d{AEMF)/dx > 0, with x 
the atomic doping content and AEmf = -^'mf ~ ^mf- This is likely because the tetragonal 
phase, which has an higher dielectric constant, better screens the electrostatic perturbation 
^doping induced by the doping^. 

In order to have a picture of the screening effect we have defined (see Fig. [5] and caption) 
the polarization charge, that is the difference in the electronic density between the pure and 
the doped system. In Fig. |5] we see that the main perturbation is induced by the oxygen 
vacancies. The screening mechanism however is very efficient and accordingly the value of 
the polarization charge is sensibly different from zero only on the atoms NN to the vacancies. 

The xc-energy in pure Zr02 instead is lower in the (M) phase and on average dAExc/ dx < 
0. To understand this point we consider the perturbations induced in the electronic part of 
the Hamiltonian of pure Zr02, that is ^/p^*"* = V^^^^g + V^^^^, where V^^^^ is the additional 
perturbation generated from the dislocation of the atoms. V^p^*** has the shape of a random 
potential in a perfect periodic system which tends to destroy the collective behavior of the 
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Figure 5. (colors online) 3D plot of the polarization charge, in the tetragonal phase, defined as: 
Ppoi = [32p(Zro.940i.97^o.06Vb.o3) + /o(0)]-32/9(Zr02): iso-surfaces for ppoi = ±po are represented, 
with the positive ones in red (lighter), the negatives in blue (darker). Oxygen atoms are the small 
spheres (red). Zirconium atoms the dark big spheres (magenta) and Yttrium atoms the light biggest 
spheres (yellow). The (green) shadowed regions highlight the regions around an oxygen vacancy 
(small black sphere in center of the region), where ppoi is different from zero. 



electron gas and thus electrons lose correlation. This is not an "on site" correlation, but 
rather a "ranged" correlation energy and so we expect that GGA can better describe this 
effect than LDA. While V^^^g is similar for the two phases, i.e. we are considering the same 
kind of doping, the V^^^^ term is greater in the (T) phase due to its higher screening (see 
note 3l|. Accordingly the xc-energy is more penalized, that is dE^tc'^ /dx < dE'ic} /dx. This 
explains dAE^c/dx < and thus we can infer that in the screening process electrons lose 
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xc-energy. 

The emerging picture is that the (M) — > (T) DIPT is a balance between different mech- 
anisms with a key role played by vacancies^, though the role of Yttrium atoms is not only 
to induce vacancies (see Fig. [Hand caption). 

III. CONCLUSIONS 

The phase transition of Zirconia induced with Yttria doping is very well known experi- 
mentally and indeed Yttria stabilized Zirconia is commonly used in many applications. Our 
results show that an ab-initio approach is able to correctly describe this mechanism, with 
a predicted phase transition at a doping concentration of ~ 7.5%, in good agreement with 
experimental data. This result confirms the opportunity to predict the effect of other less 
characterize kind of doping which could be of potential interest for new applications, as for 
example doping with magnetic materials^. 

Moreover, within density functional theory we can explore how physical properties of the 
system are influenced by doping. As an example we considered how different components 
of the energy changes and in particular the behavior of the exchange-correlation energy. 
We showed that oxygen defects play a major role in the phase transition and how the 
perturbation induced by oxygen vacancies to the system is screened in the high-k tetragonal 
phase. The same approach could be used to check the effect of doping on other quantities, 
which are not easily accessible experimentally; among others, the value of the dielectric 
constant, which has a key importance for applications in the field of micro-electronics. 

Finally these results, improving our understanding of the phase transition, can also be 
used to tune the parameters of models for the description of realistic devices based on 
Zirconia, whose dimensions are still beyond the capabilities of first-principles simulations. 
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the contribution from the lattice is dominant. The (T) phase has higher dielectric constant and 
accordingly we expect a greater distortion of the crystal, provided that the linear order behavior 
is respected to all orders. 

32 Recent experimental results seems to suggest that the (T) phase can be stabilized simply in- 
ducing oxygen vacancies, tailoring the growing conditions of the sample, without the inclusion 
of other dopants^. 
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